/*
 * main.cpp
 *
 *  Created on: 06.02.2013
 *      Author: Roma Jam
 */

#include "time.h"
#include "stdlib.h"
#include "stdio.h"
#include "inttypes.h"
#include "fft.h"

#define N0      9           // dlina registra
#define K       (1 << N0)   // dlina PSP
#define INT     2*K         // Period of search

int aPSP[INT]; // PSP
unsigned int RG = 0xffffffff;
unsigned int O5 = (1 << 5);
unsigned int O9 = (1 << 9);
unsigned int O10 = (1 << 10);
complex<double> cPSP1[INT],
                cPSP2[INT],
                cPSP3[INT],
                cPSP4[INT];
double Fd = 0;
double Fref = 0;
double dT = 2e-6;
double Ttek = 0;
double r1 = 0;
double r2 = 0;


using namespace std;
FILE *fout;

double RANDOM(){
    double r;
    do{
        r = (double)rand() / RAND_MAX ;
    } while(r == 0.0);
    return r;
}

void gen_PSP(void){
    int k;

    Ttek = (double)rand() / RAND_MAX * 12357913.0;
    for(k = 0; k<INT-1; k++){
        RG = (RG & 0xfffffffe) | (((RG & O5)>>5)^((RG & O9)>>9));
        RG <<= 1;
        aPSP[k] = (RG & O10)>>10;
        if(aPSP[k] == 0) aPSP[k] = -1;
    }
    int indx_print = 0;
    for(k = 0; k<INT-1; k++){
        r1 = RANDOM();
        r2 = RANDOM();

        cPSP1[k] = (double)aPSP[(k+100)%K] * 0.5 * exp(complex<double>(0.0,-2.0*M_PI*Fd*Ttek));
//                 + (double)aPSP[(k+400)%K] * exp(complex<double>(0.0,-2.0*M_PI*Fd*Ttek));
//                 + (complex<double>(1*sqrt(-2.0*log(r1))*cos(2.0*M_PI*r2),1*sqrt(-2.0*log(r1))*sin(2.0*M_PI*r2)));

        cPSP2[k] = (double)aPSP[k] * exp(complex<double>(0.0,-2.0*M_PI*Fref*Ttek));   // form PSP2

        cPSP3[k] = (double)aPSP[k] * exp(complex<double>(0.0,-2.0*M_PI*Fref*Ttek));         // real PSP3


        Ttek += dT;
        if(indx_print == k){
            printf(".");
            indx_print = k+21;
        }
    }

//    for(k = 0; k<K-1; k++){
//        cPSP1[(k+100)%K] += cPSP4[(k+50)%K];
//    }
}

int main(void){
    int k = 0;
    int indx_print;


    // randomize();
    srand((int)time(0));
    fout = fopen("PSP.txt","wt");
    printf("Generate PSP start");

    gen_PSP();
    fftm(cPSP1,INT,0);
    fftm(cPSP2,INT,0);
    fftm(cPSP3,INT,0);
    for(k=0;k<INT;k++){
        cPSP1[k] = complex<double>(real(cPSP1[k]),-imag(cPSP1[k]));  // made complex from PSP
        cPSP2[k] *= cPSP1[k];
        cPSP3[k] *= cPSP1[k];
    }
    fftm(cPSP2,INT,1);
    fftm(cPSP3,INT,1);

    printf("Completed!\rSaving txt");
    indx_print = 0;
    for(k=0;k<INT;k++){
        fprintf(fout,"%d %lf %lf %lf %lf\n",
                aPSP[k],
                real(cPSP2[k]), imag(cPSP2[k]),
                real(cPSP3[k]), imag(cPSP3[k])
        );
        if(indx_print == k){
            printf(".");
            indx_print = k+21;
        }
    }

    printf("'PSP.txt' has been saved.\r");
    fclose(fout);
}
